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Abstract 

A tensor network is a product of tensors associated with vertices of some graph G such that every 
edge of G represents a summation (contraction) over a matching pair of indexes. It was shown recently 
by VaHant, Cai, and Choudhary that tensor networks can be efficiently contracted on planar graphs if 
components of every tensor obey a system of quadratic equations known as matchgate identities. Such 
tensors are referred to as matchgate tensors. The present paper provides an alternative approach to 
contraction of matchgate tensor networks that easily extends to non-planar graphs. Specifically, it is 
shown that a matchgate tensor network on a graph G of genus g with n vertices can be contracted in 
time T — poly{n) + 0{m^) 2^^ where m is the minimum number of edges one has to remove from G in 
order to make it planar. Our approach makes use of anticommuting (Grassmann) variables and Gaussian 
integrals. 
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1 Introduction and summary of results 



Contraction of tensor networks is a computational problem having a variety of applications ranging from 
simulation of classical and quantum spin systems UJ 121 [21 IH |5l to computing capacity of data storage 
devices 161. Given the tremendous amount of applications it is important to identify special classes of tensor 
networks that can be contracted efficiently. For example, Markov and Shi found a linear time algorithm 
for contraction of tensor networks on trees and graphs with a bounded treewidth [IJ. An important class 
of graphs that do not fall into this category are planar graphs. Although contraction of an arbitrary tensor 
network on a planar graph is a hard problem, it has been known for a long time that the generating function 
of perfect matchings known as the matching sum can be computed efficiently on planar graphs for arbitrary 
(complex) weights using the Fisher-Kasteleyn-Temperley (FKT) method, see f71, ^Sl, IQ]. It is based on the 
observation that the matching sum can be related to Pfaffian of a weighted adjacency matrix (known as the 
Tutte matrix). The FKT method also yields an efficient algorithm for computing the partition function of 
spin models reducible to the matching sum, most notably, the Ising model on a planar graph [10]. Recently 
the FKT method has been generalized to the matching sum of non-planar graphs with a bounded genus |TT[ 
MM- 

Computing the matching sum can be regarded as a special case of a tensor network contraction. It is 
therefore desirable to characterize precisely the class of tensor networks that can be contracted efficiently 
using the FKT method. This problem has been solved by Valiant |[T4l[T5l and in the subsequent works by Cai 
and Choudhary 1.16., 17. ,181. Unfortunately, it turned out that the matching sum of planar graphs essentially 
provides the most general tensor network in this class, see lfl6l [TSl . Following |[T6ll we shall call such 
networks matchgate tensor networks, or simply matchgate networks. A surprising discovery made in fVh is 
that matchgate tensors can be characterized by a simple system of quadratic equations known as matchgate 
identities which does not make references to any graph theoretical concepts. Specifically, given a tensor T 
of rank n with complex- valued components T{x) = T'^i,x2,...,a;„ labeled by n-bit strings x G {0, 1}" one 
calls T a matchgate tensor, or simply a matchgate, if 

r(xee'')r(yee")(-lfi+-+^«-i+2/i+-+2/a-i =0 forall x.^ejO,!}". (1) 

a : Xa^Va 

Here denotes a string in which the a-th bit is 1 and all other bits are 0. The symbol © stands for a bit-wise 
XOR of binary strings. For example, a simple algebra shows that a tensor of rank n = 1, 2, 3 is a matchgate 
iff it is either even or Furthermore, an even tensor of rank 4 is a matchgate iff 

- r(oooo) r(iiii) + r(iioo) r(ooii) - r(ioio) r(oioi) + r(iooi) r(oiio) = o. (2) 

A matchgate network is a tensor network in which every tensor is a matchgate. 

The purpose of the present paper is two-fold. Firstly, we develop a formalism that allows one to per- 
form partial contractions of matchgate networks, for example, contraction of a single edge combining its 
endpoints into a single vertex. More generally, the formalism allows one to contract any connected planar 
subgraph G of the network into a single vertex u{G) by "integrating out" all internal edges of G. The 
number of parameters describing the contracted tensor assigned to u{G) is independent of the size of G. It 
depends only on the number of "external" edges connecting G to the rest of the network. This is the main 
distinction of our formalism compared to the original matchgate formalism of Valiant |[T4l . The ability to 
implement partial contractions may be useful for designing efficient parallel contraction algorithms. More 
importantly, we show that it yields a faster contraction algorithm for matchgate networks on non-planar 
graphs. 

'a tensor T is called even (odd) if T{x) = for all strings x with odd (even) Hamming weight. 
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Our formalism makes use of anticommuting (Grassmann) variables such that a tensor of rank n is repre- 
sented by a generating function of n Grassmann variables. A matchgate tensor is shown to have a Gaussian 
generating function that depends on 0{n?) parameters. The matchgate identities Eq. ([B can be described 
by a first-order differential equation making manifest their underlying symmetry. Contraction of tensors is 
equivalent to convolution of their generating functions. Contraction of matchgate tensors can be performed 
efficiently using the standard Gaussian integration technique. We use the formalism to prove that a tensor 
satisfies matchgate identities if and only if it can be represented by the matching sum on some planar graph. 
It reproduces the result obtained earlier by Cai and Choudhary ifTTlfTSl . Our approach also reveals that the 
notion of a matchgate tensor is equivalent to the one of a Gaussian operator introduced in |fT9l in the context 
of quantum computation. 

Secondly, we describe an improved algorithm for contraction of matchgate networks on non-planar 
graphs. Let S be a standard oriented closed surface of genus g, i.e., a sphere with g handles. 

Definition 1. Given a graph G = {V, E) embedded into a surface S we shall say that G is contractible if 
there exists a region D C T, with topology of a disk containing all vertices and all edges of G. A subset of 
edges M E is called a planar cut ofG if a graph Gm = (y, E\M) is contractible. 

A contraction value c{T) of a tensor network T is a complex number obtained by contracting all tensors 
of T. Our main result is as follows. 

Theorem 1. Let T be a matchgate tensor network on a graph G = iy, E) with n vertices embedded into 
a surface of genus g. Assume we are given a planar cut of G with m edges. Then the contraction value 
c{T) can be computed in time T = 0{{n + m)^) + 0{m^) 2^^. If G has a bounded vertex degree, one can 
compute c{T) in time T = 0{{n + m)^) + 0{m^) 2^^'. 

If a network has a small planar cut, m <<^n, the theorem provides a speedup for computing the matching 
sum and the partition function of the Ising model compared to the FKT method. For example, computing 
the matching sum of a graph G as above by the FKT method would require time T = 0{n^) 2^^ since the 
matching sum is expressed as a linear combination of 2^^ Pfaffians where each Pfaffian involves a matrix 
of size n X n, see l[m[T2l[T3l . and since Pfaffian of an n x n matrix can be computed in time 0(n'^), see 
Remark 2 below. In contrast to the FKT method, our algorithm is divided into two stages. At the first stage 
that requires time 0((n + m)^) one performs a partial contraction of the planar subgraph Gm determined 
by the given planar cut M, see Def . U The contraction reduces the number of edges in a network down to m 
without changing the genu^ The first stage of the algorithm yields a new network T' with a single vertex 
and m self-loops such that c(T') = c(T). At the second stage one contracts the network T' by expressing 
the contraction value c(T') as a linear combination of 2"^^ pfaffians similar to the FKT method. However 
each Pfaffian involves a matrix of size only 0{m) x 0{m). 

Remark 1: The statement of the theorem assumes that all tensors are specified by their generating functions. 
Thus a matchgate tensor of rank d can be specified by 0{d?) parameters, see Section [3] for details. The 
ordering of indexes in any tensor must be consistent with the orientation of a surface. See Section l2?T] for a 
formal definition of tensor networks. 

Remark 2: Recall that Pfaffian of an n x n antisymmetric matrix A is defined as 

Pf (A) = l if n is odd, 

I ^hj. T^aeSn sgn(o-) -4^(i),<t{2) ^<7{3),a(4) • • • ^a(n-i),<7{n) if n is even. 

where 5„ is the symmetric group and sgn((T) = ±1 is the parity of a permutation a. One can efficiently 
compute Pfaffian up to a sign using an identity Pf (A)'^ = det (A). However, in order to compute a linear 

^If the initial network represents a matchings sum, the first stage of the algoritlim would require only time 0{{n + m)^). 



4 



combination of several Pfaffians one needs to know the sign exactly. One can directly compute Pf (A) 
using the combinatorial algorithm by Mahajan et al f20l in time O(n^). Alternatively, one can use Gaussian 
elimination to find an invertible matrix U such that AU is block-diagonal with all blocks of size 2 x 2. It 
requires time 0{n^). Then Pf {A) can be computed using an identity Pf {U A U'^) = det {U) Pf {A). This 
method yields O(n^) algorithm although it is less computationally stable compared to the combinatorial 
algorithm of 1201. 

2 Some definitions and notations 

2.1 Tensor networks 

Throughout this paper a tensor of rank d is a d-dimensional complex array T in which the indexes take 
values and 1. Given a binary string of indexes x = {xiX2 ■ ■ ■ Xd) we shall denote the corresponding 
component Tr,^x2...Xi as T{x). 

A tensor network is a product of tensors whose indexes are pairwise contracted. More specifically, each 
tensor is represented by a vertex of some graph G = (V, E), where F is a set of vertices and is a set of 
edges. The graph may have self-loops and multiple edges. For every edge e £ E one defines a variable 
x{e) taking values and 1. A bit string x that assigns a particular value to every variable x(e) is called an 
index string. A set of all possible index strings will be denoted X{E). In order to define a tensor network 
on G one has to order edges incident to every vertex. We shall assume that G is specified by its incidence 
list, i.e., for every vertex u one specifies an ordered list of edges incident to u which will be denoted E{u). 
Thus E{u) = {e", . . . , e^(„)} where £ E for all j. Here d{u) = \E{u) \ is the degree of u. If a vertex u 
has one or several self-loops, we assume that every self-loop appears in the list E{u) twice (because it will 
represent contraction of two indexes). For example, a vertex with one self-loop and no other incident edges 
has degree 2. A tensor network on G is a collection of tensors T = {Tu}ueV labeled by vertices of G such 
that a tensor has rank d{u). A contraction value of a network T is defined as 

^('^)= E n Tu{x{e\)...x{e\^~^)). (3) 
xgx{e) uev 

Thus the contraction value can be computed by taking a tensor product of all tensors {T^} and then contract- 
ing those pairs of indexes that correspond to the same edge of the graph. By definition, c(T) is a complex 
number (tensor of rank 0). 

It will be implicitly assumed throughout this paper that a tensor network is defined on a graph G embedded 
into a closed oriented surface E. We require that the order of edges incident to any vertex u must agree with 
the order in which the edges appear if one circumnavigates u counterclockwise. Thus the order on any set 
E{u) is completely specified by the choice of the first edge e" G E{u). If the surface S has genus g we 
shall say that G has genus g (it may or may not be the minimal genus for which the embedding of G into S 
is possible). 

2.2 Anticommuting variables 

In this section we introduce notations pertaining to the Grassmann algebra and anticommuting variables 
(see the textbook ET\ for more details). Consider a set of formal variables 9 = {9i, . . . ,9n) subject to 
multiplication rules 

9l = 0, 9a9b + 9b9a = for all a,b. (4) 
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The Grassmann algebra Q{9) is the algebra of complex polynomials in variables 9i, . . . ,9^ factorized over 
the ideal generated by Eq. (01). Equivalently, Q{9) is the exterior algebra of the vector space C", where each 
variable 9a is regarded as a basis vector of C". More generally, the variables 6a may be labeled by elements 
of an arbitrary finite set X (in our case the variables will be associated with edges or vertices of a graph). A 
linear basis of Q{6) is spanned by 2" monomials in variables Oa- Namely, for any subset M C {1, . . . , n} 
define a normally ordered monomial 

9{M) =\{ea (5) 

aeM 

where the indexes increase from the left to the right. If the variables are labeled by elements of some set X, 
one can define the normally ordered monomials 9{M), M C X by choosing some order on X. Let us agree 
that ^(0) = /. Then an arbitrary element f £ Q[9) can be written as 

/= Yl f(M)eiM), /(M)€C. (6) 

MC{l,...,n} 

We shall use notations / and f{9) interchangeably meaning that / can be regarded as a function of anticom- 
muting variables 9 = (^i, . . . , Accordingly, elements of the Grassmann algebra will be referred to as 
functions. In particular, / is regarded as a constant function. A function f{9) is called even (odd) if it is a 
linear combination of monomials 6{M) with even (odd) degree. Even functions span the central subalgebra 

of g{9). 

We shall often consider several species of Grassmann variables, for example, 6 = (6'i, . . . , 6'„) and 
r/ = {rji, . . . ,r]k). It is always understood that different variables anticommute. For example, a function 
f{9, T]) must be regarded as an element of the Grassmann algebra Q{9, rj), that is, a linear combination of 
monomials in ^i, . . . , and 771, ... , r]k. 

A partial derivative over a variable 9a is a linear map da : Q{9) ^ Q{0) defined by requirement da-I = Q 
and the Leibniz rule 

da-{ebf) = 5a,bf -eb{da- f). 

More exphcitly, given any function / € G{9), represent it as f{0) = fo + 9a fi, where /o, /i G G{9) do not 
depend on 9a- Then da f = fi- It follows that da • 9a = I, da9h = —9bda, dadb = —dbda for a ^ b and 
dl = 0. 

A linear change of variables 9a = Ylb=i ^a,b ^fe with invertible matrix U induces an automorphism of the 
algebra G{9) such that f{6) —>■ f{9). The corresponding transformation of partial derivatives is 

n 

da = Y.^U-\,adb. (7) 
6=1 



2.3 Gaussian integrals 

Let 9 = (6*1, . . . , 0„) be a set of Grassmann variables. An integral over a variable Oa denoted by J d9a is 
a linear map from Q{9i, . . . ,9n) to Q{9i, . . . ,9a, ■ ■ ■ ,9n), where 9 a means that the variable 9 a is omitted. To 
define an integral J d9a f{9), represent the function / as / = /0+6'a fi, where /o, /i G G{9i, . . . ,9a, ■ ■ ■ , 9^). 
Then J d9af{9) = fi. Thus one can compute the integral f d9af{9) by first computing the derivative 
da ■ f{9) and then excluding the variable 9a from the list of variables of /. 

Given an ordered set of Grassmann variables 9 = (^1, . . . , 0„) we shall use a shorthand notation 

D9= f d9n--- f d92 [ d9i. 
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Thus f DO can be regarded as a linear functional on G{0), or as a linear map from Q{9,r]) to G{r]), and so 
on. The action of / DO on the normally ordered monomials is as follows 

De9{M) = ll M = {l,2,...,n}, 
^ ^ 10 otherwise. 



Similarly, if one regards J D9 as a linear map from Q{9, rj) to G{i]) then 

DOO{M)ri{K) = 



I 



r^{K) if M = {l,2,...,n}, 
otherwise. 



Although this definition assumes that both variables 9, rj have a normal ordering, the integral / DO depends 
only on the ordering of 0. 

One can easily check that integrals over different variables anticommute, J dOa J dOf, = — J dOf, J dOa 
for a ^ b. More generally, if 6* = (6*1, ... , On) and rj = (r/i, . . . , rjk) then 



DO j Di] = (-l)"'^ j J DO. (9) 

Under a linear change of variables 9a = Ylb=i ^a,b Vb the integral transforms as 

j D0 = det ([/) j Drj. (10) 

In the rest of the section we consider two species of Grassmann variables = {Oi, . . . ,9n) and rj = 
. . . , r/fc). Given an antisymmetric n x n matrix A and any n x k matrix B, define quadratic forms 

n n k 

9'^ A9 = Aa,b9a9b, O^Bri = Y,Y.^-^bOaVb- 

a,b=l a=l b=l 

Gaussian integrals over Grassmann variables are defined as follows. 

I{A) = j D9 exp Q 0^ A 9^ and I{A, B) = j D9 exp (^9'^ AO + O'^ B . (11) 

Thus I{A) is just a complex number while I{A, B) is an element of G{ti). Below we present the standard 
formulas for the Gaussian integrals. Firstly, 

/(^)=Pf(^). (12) 

Secondly, if A is an invertible matrix then 

I{A, B) = Pf {A) exp Q rf B'^A'^B . (13) 

Assume now that A has rank m for some eveiH integer < m < n. Choose any invertible matrix U such 
that AU has zero columns m + 1, . . . , n. (This is equivalent to finding a basis of C" such that the last n — m 
basis vectors belong to the zero subspace of ^4.) Then 



U'^ AU 



All 




^Note that antisymmetric matrices always have even rank. 
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for some invertible ni x m matrix An. Introduce also matrices Bi, B2 of size m x k and {n — m) x k 
respectively such that 



B 



Bo 



Performing a change of variables 9 = U6 in Eq. ([TTI ) and introducing variables r = (ti, . . . ,rm) and 
/i = (/XI, . . . , fin-m) such that 9 = (r, /i) one gets 

I{A, B) = det {U) j Dt exp Q An t + t'^ Bi ^ exp (/i^ ^2 r/) . 

Here we have taken into account Eqs. (I9I10I ). Applying Eq. (fT3l) to the first integral one gets 

/(AS) = Pf(^ii) det(C/) expQ?7^Sf(An)-'Bir?) y" I?/x exp (/x^ r?) . (14) 

One can easily check that / -D/U exp (/i^ i?2 =0 if the rank of B2 is smaller than the number of variables 
in /i, that is, n — m. Since B2 has only A; columns we conclude that 

I{A, B) = unless m > n — k. 

Therefore in the non-trivial case I{A, B) ^ the matrices Bf {Au)~^ Bi and B2 specifying I{A, B) have 
size k X k and k' x k for some k' < k. It means that I{A, B) can be specified by 0{k'^) bits. One can 
compute I{A,B) in time 0(n^ + v?k). Indeed, one can use Gaussian elimination to find U, compute 
det {U) and Pf (An) in time O(n^). The matrix A'^\ can be computed in time O(n^). Computing the 
matrices i?2 requires time 0{'n?k). 

The formula Eq. (fT4l ) will be our main tool for contraction of matchgate tensor networks. 



3 Matchgate tensors 

3.1 Basic properties of matchgate tensors 

Although the definition of a matchgate tensor in terms of the matchgate identities Eq. ([T|l is very simple, it 
is neither very insightful nor very useful. Two equivalent but more operational definitions will be given in 
Sections [331 13.41 Here we list some basic properties of matchgate tensors that can be derived directly from 
Eq. ([U. In particular, following the approach of \XT\, we prove that a matchgate tensor of rank n can be 
specified by a mean vector z E {0, 1}" and a covariance matrix A of size n x n. 

Proposition 1. Let T be a matchgate tensor of rank n. For any z € {0, 1}" a tensor T' with components 
T'{x) = T{x © z) is a matchgate tensor 

Proof. Indeed, make a change of variables x^x(Bz,y^y®zm the matchgate identities □ 

Let r be a non-zero matchgate tensor of rank n. Choose any string z such that T{z) ^ and define a 
new tensor T' with components 

r'(x) = ^|^, xG{o,ir, 

such that T' is a matchgate and T'(0") = 1. Introduce an antisymmetric n x n matrix A such that 

r'(e" © a^) if a<b, 
Aa,b = { -T'ie" © a^) if a > b, 
if a = b. 
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Proposition 2. For any x € {0, 1}" 

Vi{A{x)) if X has even weight 
\ if X has odd weight ' 

where A{x) is a matrix obtained from A by removing all rows and columns a such that Xa = 0. 

Proof. Let us prove the proposition by induction in the weight of x. Choosing x = 0" and y = e"' in the 
matchgate identities Eq. ^ one gets r'(e") = for all a. Similarly, choosing x = and y = e°- with 
a <h one gets T'{e°- © e^) = ^a,f) = Pf {A{e"- © e*)). Thus the proposition is true for |x| = 1, 2. Assume 
it is true for all strings x of weight < k. For any string x of weight k + 1 and any a such that = apply 
the matchgate identities Eq. ([ij with x and y = e". After simple algebra one gets 

6-1 

r'(x©e'^)= ^a,feT'(x©e^)(-ir("'''), r?(a, 6) = J] x,-. 

h:xi,=l j=a 

Noting that x © has weight k and applying the induction hypothesis one gets 

r'(x©e")= ^a,6 Pf(^(x©e^)) 

b:xt = l 

for even k and T'{x ffi e") = for odd k. Thus T'(y) = for all odd strings of weight k + 2. Furthermore, 
let non-zero bits of x © be located at positions ji < j2 < • • • < jk- Note that the sign of Aa^b (—1 
coincides with the parity of a permutation that orders elements in a set [a, b,ji,j2, . . . ,jk]- Therefore, by 
definition of Pfaffian one gets T'{x © e") = Pf {A{x © e"")). □ 

Thus one can regard the vector z and the matrix A above as analogues of a mean vector and a covariance 
matrix for Gaussian states of fermionic modes, see for instance [19|. Although Proposition [2] provides 
a concise description of a matchgate tensor, it is not very convenient for contracting matchgate networks 
because the mean vector z and the covariance matrix A are not uniquely defined. 

Corollary 1. Any matchgate tensor is either even or odd. 

Proof. Indeed, the proposition above implies that if a matchgate tensor T has even (odd) mean vector it is 
an even (odd) tensor. □ 



3.2 Describing a tensor by a generating function 

Let 6 = [01, ... , On) be an ordered set of n Grassmann variables. For any tensor T of rank n define a 
generating function T G Q(9) according to 

T{e)= Y T{x)e{x). 

xe{o,i}" 

Here 6{x) = 6^^ ■ ■ ■ is the normally ordered monomial corresponding to the subset of indexes x = {a : 
Xa = 1}. Let us introduce a linear differential operator A acting on the tensor product of two Grassmann 
algebras Q{0) ® Q{0) such that 

n 
a=l 



9 



Lemma 1. A tensor T of rank n is a matchgate iff 



A-T(^T = 0. (16) 
Proof. For any strings x,y € {0, 1}" one has the following identity: 

Expanding both factors T in Eq. ([T6l ) in the monomials 9{x), 9{y), using the above identity, and performing 
a change of variable X — > X Cq, and y — > y Ca for every a one gets a linear combination of monomials 
9{x) (g) 9{y) with the coefficients given by the right hand side of Eq. Therefore Eq. ([T6l ) is equivalent to 
Eq. ©. □ 

Lemma [T] provides an alternative definition of a matchgate tensor which is much more useful than the 
original definition Eq. ([T]). For example, it is shown below that the operator A has a lot of symmetries which 
can be translated into a group of transformations preserving the subset of matchagate tensors. 

Lemma 2. The operator A is invariant under linear reversible changes of variables 9a = ^5=1 Ua,b 9b- 
Proof. Indeed, let da be the partial derivative over 9a. Using Eq. (|7]) one gets 

n n n 

^ea®da + da®9a= ^^^b {U-^)c,a {h ®dc + dc®h)=Y^ {h + h)- 

a=l a,b,c=l b 

□ 

Lemmas 1 1 121 imply that Unear reversible change of variables T{9) T{9), where 9a = Y^^=i Ua,bGb 
map matchgates to matchgates. 

Corollary 2. Let T be a matchgate tensor of rank n. Then a tensor T' defined by any of the following 

transformations is also matchgate. 

(Cyclic shift): T'{xi,X2, ■■■ ,Xn) = T{x2, • • • , x„, xi), 

(Reflection): T'{xi,X2, ■ ■ ■ ,Xn) = T{xn, ■ ■ ■ ,X2,xi), 

(Phase shift): T'{x) = (-l)^'^ T{x), where z G {0, 1}". 



Proof. Let e = if T is an even tensor and e = 1 if T is an odd tensor, see Corollary [T] The transformations 
listed above are generated by the following linear changes of variables: 

Phase shift : 9a (-1)^" 9a, a = 1, . . . , n. 
Cyclic shift : 9a 9a-i a = 2, . . . , n, and ^ (-1)^^^ ^n- 
Reflection : 9a^i9n-a- 

Indeed, let 9{x) be the normally ordered monomial where x = (xi, X2, • • • , x„). Let x' = {x2, ■ ■ ■ , Xn, xi) 
for the cyclic shift and x' = . . . ,X2,xi) for the reflection. Then the linear changes of variables stated 
above map 9{x) to (—1)^'^ 9{x) for the phase shift, to 9{x') for the cyclic shift, and to i'^9{x') for the 
reflection. Therefore, in all three cases T' is a matchgate tensor. □ 
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3.3 Matchgate tensors have Gaussian generating function 

A memory size required to store a tensor of rank n typically grows exponentially with n. However the 
following theorem shows that for matchgate tensors the situation is much better. 

Theorem 2. A tensor T of rank n is a matchgate iff there exist an integer < k < n, complex matrices A, 
B of size n X n and k x n respectively, and a complex number C such that T has generating function 

T{9) = C exp Q 0^ A j Dfi exp (^^ B 6) , (17) 

where jj, = (/xi, . . . is a set of k Grassmann variables. Furthermore, one can always choose the 
matrices A and B such that = —A and BA = 0. 

Thus the triple {A, B, C) provides a concise description of a matchgate tensor that requires a memory 
size only O(n^). In addition, it will be shown that contraction of matchgate tensors can be efficiently 
implemented using the representation Eq. ( fTT] ) and the Gaussian integral formulas of Section [231 We shall 
refer to the generating function Eq. ([TT] ) as a canonical generating function for a matchgate tensor T. 

Corollary 3. For any matrices A and B the Gaussian integral I{A, B) defined in Eq. ([77]) is a matchgate. 

Proof. Indeed, use Eq. ([T4l) and Theorem [2] □ 

In the rest of the section we shall prove Theorem |2l 

Proof of Theorem^ Let us first verify that the tensor defined in Eq. (fT/l) is a matchgate, i.e., A • T ® T = 0, 
see Lemma[I] Without loss of generality A is an antisymmetric matrix and C = 1. Write T as 

T = T2Ti, where Ta = exp Q 6*^ ^ 6*^ , Ti = j D fi exp {fi'^ B 9) . 

Noting that T2 is an even function and da 9{x) = da ■ 9{x) + 9{x) da for any even string x one concludes 
that 

A • T ® T = (A • Ta ^ Ta) Ti Ti + Ta ® T2 (A • Ti Ti) . (18) 

Therefore it suffices to prove that A • Ta Ta = and A • Ti (g) Ti = 0. The first identity follows from 
da-T2 = ELi ^a,b 9b T2 and A^ = -A which imphes 

n 

A • Ta ® Ta = ^ Aa,b {9a ®9b + 9b® 9a) Ta O Ta = 0. 

a, 6=1 

To prove the second identity consider the singular value decomposition B = L?- BR, where L G SU{k) 
and R € SU{n) are unitary operators, while B is ak x n matrix with all non-zero elements located on the 
main diagonal, B = diag(i?i, . . . , Bj^). Introducing new variables 9 = R9 and fi = L fi one gets 

Ti = J Dfiexp (^^Bajla9~}j = Bi ■ ■ ■ bJi ■ ■ ■ 9k. 

Here we have used identity / Dfi = det (L) / Dfi = f Dfi, see Eq. ([TOl l. Since A is invariant under linear 
reversible changes of variables, see Lemma |2l and since A ■ 9{x) (S> 9{x) = for any monomial 9{x) one 
gets A • Ti (g) Ti = 0. We proved that A • T ® T = 0, that is, T is a matchgate tensor. 
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Let us now show that any matchgate tensor T of rank n can be written as in Eq. (fTTl ). Define a linear 
subspace Z C C" such that 

n 

Z = {eGC" : Y.UOaT = 0}. 

a=l 

Let dim {Z) = k. Make a change of variables rj = U 6 where U is any invertible matrix such that the last k 
rows of U span Z. Then rya T = for all a = n — /c + 1, n. It follows that T can be represented as 

T = rjn-k+i ■■■nnS (19) 
for some function S = S{r]) that depends only on variables 771, ... , rjn-k- Equivalently, 

S = dn - ■ ■ dn-k+1 ■ T, 

where the partial derivatives are taken with respect to the variables r]. Since A is invariant under reversible 
linear changes of variables, see Lemma|2j and since A da (E> da = da (S> da A, we get 

n—k 

A- S(^S = ^r,aS®da- S + da-S®r]aS = 0. (20) 

a=l 

By definition of the subspace Z the functions rjiS, . . . ,r]n-kS are linearly independent. Therefore there 
exist linear functional Fa : Q{rj) ^ C, a = 1, . . . ,n — k, such that Fa{r)bS) = 6a,b- Applying Fa to the 
first factor in Eq. (l20l ) we get 

n—k 

da-S = Y,Ma,bVbS, Where Ma,b = -Fa{db ■ S) e C, (21) 
b=i 

for all a = I, . . . ,n — k. Let kmin the lowest degree of monomials in S. Let us show that kmin = 0, that 
is, S{r]) contains / with a non-zero coefficient. Indeed, let Smin be a function obtained from S by retaining 
only monomials of degree kmin- Since any monomial in the r.h.s. of Eq. (|2TI) has degree at least kmin + 1> 
we conclude that da ■ Smin = for all a. It means that Smin = C / for some complex number C 7^ and 

thus kmin = 0. 

Applying the partial derivative d}, to Eq. (|2TI ) we get 5 = C^^{di, da ■ S)\^^q, where the substitution 
r] = means that the term proportional to the identity is taken. Since the partial derivatives over different 
variables anticommute, M is an antisymmetric matrix. 

Using Gaussian elimination any antisymmetric matrix M can be brought into a block-diagonal form with 
2x2 blocks on the diagonal by a transformation M — > M' = X W, where W is an invertible matrix 
(in fact, one can always choose unitary W, see ll23l ). Since our change of variables t] = U6 allows arbitrary 
transformations in the subspace of r/i, . . . , rjn-k we can assume that M is already bock-diagonal. 



a=l 

where only non-zero blocks are represented, so that 2m < n — k. 
Applying Eq. (|2TI ) for o = 1, 2 we get 

di-S = XimS, d2-S = -XimS. (22) 

Note that S can be written as 

S = "^{axm + PxV2)v{x) + ^{lyl + ^ymV2)r]{y), (23) 



12 



where the sums over x and y run over all odd and even monomials in 7/3, ... , i^n-k respectively. Substituting 
Eq. (I23] ) into Eq. (|22)) one gets ax = fix = ^ and 5x = Ai72:, that is 

5 = (/ + Ai??ir/2)5', 

where S" depends only on variables rjs, . . . ,r]n-k- Repeating this argument inductively, we arrive to the 
representation 

m /I \ 

S = Cl[{I + XaV2a-lV2a) = C exp ( -7?^ M 7? J . 
a=l ^ ^ 

Here we extended the matrix M such that its last k columns and rows are zero. Combining it with Eq. ([T9b 
one gets 

T = C rin^k+i ■■■r]n exp Qr;'^ M r]^ = C exp (^rf M i]^ j Dfi exp (^/x^ S r/ 
where /i is a vector of /c Grassmann variables and B is a.k x n matrix with 0,1 entries such that 

k 

f/' Br] = '^fia Vn-k+a- 
a=l 

Recalling that r] = U 9, we conclude that T has a representation Eq. ([17] ) with A = M U and B = BU. 
As a byproduct we also proved that the matrices A, B in Eq. (fTTl) can always be chosen such that BA = 
since BA = B MU and all non-zero entries of B are in the last k rows. □ 



3.4 Graph theoretic definition of matchgate tensors 

Let G = {V,E,W)he. nn arbitrary weighted graph with a set of vertices V , set of edges E and a weight 
function W that assigns a complex weight W{e) to every edge e £ E. 

Definition 2. Let G = {V, E) be a graph and S" C V be a subset of vertices. A subset of edges M E 
is called an S -imperfect matching iff every vertex from S has no incident edges from M while every vertex 
from V\S has exactly one incident edge from M. A set of all S -imperfect matchings in a graph G will be 
denoted M{G, S). 

Note that a perfect matching corresponds to an 0-imperfect matching. Occasionally we shall denote a set 
of perfect matching by M.{G) = M.{G, 0). For any subset of vertices 5 C F define a matching sum 

PerfMatch(G, 5) = ^ H 

M£M(G,S) eeM 

(A matching sum can be identified with a planar matchgate of fT5l|.) In this section we outline an isomor- 
phism between matchgate tensors and matching sums of planar graphs discovered earlier in ifTSl . For the 
sake of completeness we provide a proof of this result below. Although the main idea of the proof is the 
same as in lITSi some technical details are different. In particular, we use much simpler crossing gadget. 

Specifically, we shall consider planar weighted graphs G = {V, E, W) embedded into a disk such that 
some subset of n external vertices Vext ^ V belongs to the boundary of disk while all other internal 
vertices V\Vext belong to the interior of D. Let Vext = {^^i, . . . , be an ordered list of external vertices 
corresponding to circumnavigating anticlockwise the boundary of the disk. Then any binary string x G 
{0, 1}" can be identified with a subset x C Vext that includes all external vertices uj such that Xj = 1. Now 
we are ready to state the main result of this section. 
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Figure 1: Left: a complete graph Cq embedded into a disk. Right: a perfect matching on Cg with two 
self -intersections . 

Theorem 3. For any matchgate tensor T of rank n there exists a planar weighted graph G = (Vi W) 
with 0{n?) vertices, 0{n?) edges and a subset of n vertices Vext ^ ^ ^iich that 

r(x) = PerfMatch(G, x) for all x^Vext- (25) 

Furthermore, suppose T is specified by its generating function, T = C exp [^O^^ A9^ f Dfi exp (^fi^ B Oy 
Then the graph G can be constructed in time 0{n'^) and the weights W{e) are linear functionals of A, B, 
and C. 

The key step in proving the theorem is to show that Pfaffian of any n x n antisymmetric matrix can be 
expressed as a matching sum on some planar graph with O(n^) vertices. This step can be regarded as a 
reversal of the FKT method that allows one to represent the matching sum of a planar graph as Pfaffian of 
the Tutte matrix. 

Lemma 3. For any complex antisymmetric matrix A of size n x n there exists a planar weighted graph 
G = {V, E, W) with 0{n?) vertices, 0{'n?) edges such that the weights W{e) are linear functionals of A 
and 

Pf {A) = PerfMatch(G, 0). (26) 
The graph G can be constructed in time 0(in?). 

Remark: It should be emphasized that we regard both sides of Eq. (l26l ) as polynomial functions of matrix 
elements of A, and the lemma states that the two polynomials coincide. However, even if one treats both 
sides of Eq. (1261 ) just as complex numbers, the statement of the lemma is still non-trivial, since one can not 
compute Pf {A) in time O(n^) and thus one has to construct the graph G without access to the value of 
Pf(A). 

Proof. Let us assume that n is even (otherwise the statement is trivial). Let D be a disk with n marked 
points f 1 , . . . , u„ on the boundary such that their order corresponds to anticlockwise circumnavigating the 
boundary of D. Let C„ be the complete graph with vertices vi, . . . ,Vn embedded into D. We assume that 
the embedding is chosen such that all edges of C„ lie inside the disk and there are only double edge crossing 
points, see Fig. 1. Let A^(C„) be a set of perfect matchings on C„. For any perfect matching M € A1(C„) 
let Nc{M) be the number of self-intersections in M, i.e., the number of edge crossing points in the planar 
embedding of C„ in which both crossing edges ai^e occupied by M. For example, given a planar embedding 
of Ce shown on Fig. 1, a perfect matching M = (1, 3), (2, 5), (4, 6) has two self-intersections. We claim 
that 

pf(^)= i-^f^^'"^ n ^27) 

M&M{Cn) {u,v)eM,u<v 
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Indeed, by definition of Pfaffian 

(28) 

where the sum is over all permutations of n elements a such that a{2j — 1) < a{2j) for all j and a{l) < 
(t(3) < . . . < cr(n — 1). Clearly, there exists a one-to-one correspondence between such permutations 
and perfect matchings in C„. If M is the perfect matching corresponding to the identity permutation, 
M = (1, 2), . . . , (n — 1, n), one has Nc{M) = and the signs in Eqs. (I27I28I ) coincide. Furthermore, 
changing M by any transposition j ^ j + 1 either does not change M or changes the parity of Nc{M), so 
the signs in Eqs. (I27I28I ) coincide for all perfect matchings. 

In order to represent the sum over perfect matchings in Eq. (|27] ) as a sum over perfect matchings in a 
planar graph we shall replace each edge crossing point of C„ by a crossing gadget, see Fig. 2. A crossing 
gadget is a planar simulator for an edge crossing point. It allows one to establish a correspondence between 
subsets of edges in the non-planar graph and subsets of edges in a planar graph. In addition, a crossing 
gadget will take care of the extra sigr0 factor in Eq. (|27] ). 

Crossing gadget. Consider a weighted graph Gcross shown on Fig. 2. It has 6 vertices and 7 edges. The 
edge (5,6) carries weight —1 and all other edges carry weight +1. We fix the embedding of Gcross into 
a disk such that Gcross has four external vertices {1, 2, 3, 4} on the boundary of the disk. One can easily 
check that the matching sum of Gcross satisfies the following identities: 

PerfMatch(Gcross,0) = 1, 

PerfMatch(Gcross, {1, 3}) = PerfMatch(Gero.s, {2, 4}) = 1, 

PerfMatch(Gc^o,„{l,2,3,4}) = -1, 

PerfMatch (Gcross, {1, 2}) = PerfMatch(Gcross, {3, 4}) = 0, 

PerfMatch (Gcross, {1, 4}) = PerfMatch(Gcross, {2, 3}) = 0. 

These identities are illustrated in Fig. 3. In addition, PerfMatch (Gcross, S") = whenever 15*1 is odd. 
Thus the four boundary conditions for which the matching sum is non-zero represents the four possible 
configurations (empty/occupied) of a pair of crossing edges if they were attached to the vertices {1, 2, 3, 4}. 
For every edge crossing point of G„ one has to cut out a small disk centered at the crossing point and replace 
the interior of the disk by the gadget Gcross such that the four vertices {1, 2, 3, 4} are attached to the four 
external edges, see Fig. 2. Let G„ be the resulting graph. By construction, G„ is planar. It remains to assign 
weights to edges of G„ such that 

PerfMatch(G„, 0) = Pf (A). (29) 

Any edge of G^ falls into one of the four categories: (i) edge of G„; (ii) a section of some edge of Gn 
between two crossing gadgets; (iii) a section of some edge of Gn between a vertex of G„ and some crossing 
gadget; (iv) an edge that belongs to some crossing gadget. Note that the edges of type (iv) have been already 
assigned a weight, whereas any edge of type (i),(ii), and (iii) has a unique ancestor edge e = {u, v) in G„. 
Let us agree that for every edge e = {u, v), u < v of G„ we choose one of its descendants e in G^ and 
assign e the weight A^^v, while all other descendants of e are assigned the weight 1. Since all descendants of 
e appear or do not appear in any perfect matching M G 7W(G„) simultaneously, we arrive to Eq. (l29l ). that 
is, Cn is the desired graph G. It remains to count the number of vertices in G„. There are 0{v?) crossing 
gadgets each having 0(1) vertices. Thus G„ has O(n^) vertices. Since G„ is a planar graph it has O(n^) 
edges, see ll24l . □ 



""One can gain some intuition about the extra sign factor in Eq. l l27t if one thinks about the set of edges occupied by a perfect 
matching y as a family of "world lines" of fermionic particles. The contribution from y to Pf [A) can be thought of as a quantum 
amplitude assigned to this family of world lines. Whenever two particles are exchanged the amplitude acquires an extra factor —1. 
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Figure 2: Each edge crossing point in the planar embedding of the complete graph C„ is replaced by the 
crossing gadget G cross- Edges labeled by it carry a weight ±1. 




Figure 3: Matching sums of the graph G cross corresponding to various boundary conditions. 



Let Cn be a planar graph constructed above. Consider a matching sum PerfMatch(C'„, 5) for some 
subset C {1, 2, . . . , n} of vertices lying on the boundary of the disk. By repeating the arguments used in 
the proof of Lemma |3] one concludes that 



PerfMatch(C„,5) =Pf (^[5]) for all 5 C {1, 2, . . . , n}, 



(30) 



where A[S] is a matrix obtained from A by removing all rows and columns a € 5". Theorem|3]follows from 
Eq. (|30l ) and the following simple observation. 



Lemma 4. Let T be a matchgate tensor of rank n with a parity e(T) specified by its generating function 

T = Cexp(^^9'^Fe^ j Dfiexp{f/^G9). (31) 



Then 

T{x) = Ce(T) Pf (^A{x l'^"")) for all x G {0, 1}'", where A 
The matrix A has size k x k with n < k < 2n. 



F -G^ 
G 



(32) 



Remark: As usual, A{y) denotes a matrix obtained from A by removing all columns and rows a such that 
t/a = 0. We assume that e(r) = 1 (e(T) = — 1) for even (odd) tensors. 

Proof. Theorem |2] asserts that T always has a generating function Eq. (|3TI ) where G has size m x n for some 
m < n. Thus k = n + m < 2n. Introducing a set of k Grassmann variables t] = {6i, . . . ,9^ Hi, • • • , Hm) 
one can rewrite T as 
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Expanding the exponent one gets 



Note that 




Taking into account that m is even (odd) for even (odd) tensors and so is zi + • • • + z„ we conclude that 

T{9) = Ce{T) Pf(^(xl''""))e(x), (33) 

a;e{0,l}" 

that is r(x) = Ce{T) Pf (yl(xl'=-")). □ 

Proof of Theorem\3l Let A be the k x k matrix constructed in Lemma |4] and Cj. be the weighted planar 
graph constructed in Lemma[3]such that Eq. ( [30l ) holds for all S C {1, 2, . . . , k}. Therefore, 

T{x) = Ce{T) PerfMatch(C'fc, xO'=~") for all x G {0, 1}" (34) 

where x is obtained from x by flipping every bit. In order to transform Eq. (l34l ) into Eq. (1251 ) one can 
incorporate the factor Ce{T) into the matching sum by introducing an extra edge with a weight Ce{T) and 
adding one extra edge with weight 1 to every vertex 1, 2, . . . , n of the graph Ck in order to flip bits of x. □ 

Although it is not necessary, let us mention that the reverse of Theorem [3] is also true, namely, a tensor 
T defined by Eq. (l25l) is always a matchgate. The easiest way to prove it is to represent the matching sum 
PerfMatch(G, x) in Eq. (l25l) as a contraction of an open matchgate tensor network, see Section 14.31 in 
which every tensor has a linear generating function (thus simulating the perfect matching condition). Then 
one can use Corollary |4] to prove that T is a matchgate. 



4 Contraction of matchgate tensor networks 
4.1 Edge contractions 

Consider a tensor network T defined on a graph G = {V, E) embedded to a surface S. Suppose one can 
find a region D C S with topology of a disk such that D contains exactly two vertices u^v and several 
edges connecting u and v as shown on Fig. 4. We shall define a new tensor network T' such that: (i) T' 
coincides with T outside D; (ii) T' contains only one vertex inside D; (iii) conti^action values of T and T' 
are the same. The operation of replacing T by T' will be referred to as an edge contraction. The new vertex 
obtained by contracting all edges connecting u and v inside D will be denoted u ★ 

Suppose there are b edges connecting u and v that lie inside the disk. Applying, if necessary, a cyclic 
shift of components to the tensors Tu and/or Ty we can assume that these edges correspond to the last b 
components of the tensor Tu and the first b components of T^, see Fig. 4. Note that if the tensors under 
consideration are matchgates, the tensors obtained after the cyclic shift are also matchgates, see Corollary |2l 
In the new network T' a pair of vertices u, v is replaced by a single vertex u-k v with degree d{u -k v) = 
d{u) + d{v) — 2b. We define a new tensor Tui,v as 

Tui,v{x,y) = Y Tu{x,Zb,Zb^i,...,zi)T,a{zi,...,Zb-i,Zb,y), (35) 

zi,...,Zf,=0,l 
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D 



Figure 4: The ordering of edges before and after contraction of u and v. 




Figure 5: Contraction of self-loops can be reduced to edge contraction by adding dummy vertices. 

where x and y can be arbitrary bit strings of length d{u) — b and d{v) — b respectively. By definition of the 
contraction value, c(T) = c{T'). 

We shall also define a self-loop contraction as a special case of edge contraction. Namely, suppose one 
can find a region L) C S with topology of a disk such that D contains exactly one vertex u and several 
self-loops as shown on Fig. 5. We shall define a new tensor network T' such that: (i) T' coincides with 
T outside D; (ii) T' contains one vertex without self-loops inside D; (iii) contraction values of T and T' 
are the same. The operation of replacing T by T' will be referred to as a self-loop contraction. To define 
this operation, choose the most inner self-loop 7 G E{u) introduce a dummy vertex v near the median of 7 
and assign a tensor Tt,(xi, X2) = 5x^^x2 to this vertex. Clearly it does not change a contraction value of a 
network. Secondly, apply the edge contraction described above to the two edges connecting u and v. This 
reduces the number of self-loops by one. Repeat these two steps until all self-loops inside D are contracted. 

It should be mentioned that self-loops j £ E{u) can be identified with elements of the fundamental group 
[7] G 7ri(S,n) of the surface S with a base point u. We do not allow to contract self-loops representing 
non-trivial homotopy classes (because it cannot be done efficiently for matchgate tensor networks). 

4.2 Edge contraction as a convolution of generating functions 

Let T = {Tu}ueV be a tensor network considered in the previous section. In order to describe each tensor 
Tu by a generating function Tu{6) we shall introduce Grassmann variables 6u,i, ■ • • , Ou,d{u) associated with 
the edges e", . . . , e^^^^ G E{u) incident to u such that 

T^{e)= Yl T{x){eu,ir{eu,2r---{eu,nr", n^d{u). (36) 

a;G{0,l}" 

Similarly one can describe the contracted tensor Tui,v in Eq. (|35] | by a generating function 

x£{0,l}P J/G{0,1}« 

where p = d{u) — b and q = d{v) — h. The goal of this section is to represent the function Tui,v{G) as an 
integral of Tu{9)Ty{9) in which all variables associated with the edges to be contracted are integrated out. 
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Let E{u, v) be a set of edges connecting u and v. For any edge e G E{u, v) such that e is labeled as 

G -E(n) and as e\ G -E'(t') denote 

^(e) = ^«,,^t,,fc, / De{e)= I de.^u I dOuj, and [ De{e) = J] / D9{e). 

Note that these definitions make sense only (n, t>) is regai^ded as an ordered pair of vertices. Also note that 
the integrals J D9{e) over different edges commute, see Eq. and thus one can take the integrals in an 
arbitrary order. 

Lemma 5. Suppose the edges connecting u and v are ordered as shown on Fig. 4, i.e., these are the last b 
edges incident to u and the first b edges incident to v. Then 



Tu*v= [ D9{e)Tunexpl 0(e) 



(38) 



Proof. By linearity it is enough to prove Eq. (1381 ) for the case when and T„ are monomials in the Grass- 
mann variables, i.e., 

Tu = {0^,1^ ■ ■ ■ {Ou,pr''{0u,p+ir'' ■ ■ ■ {eu,p+br'^, % = • • • {0v,by'{0v,b+iY' ■ ■ ■ {0,,,+^?", 

where p = d{u) — b and q = d(v) — b. By expanding the exponent one gets a sum of all possible monomials 
in which the two variables associated with any edge e G E{u, v) are either both present or both absent. 
Therefore the integral in Eq. (1381 ) is zero unless Zj = for all j = 1, . . . , b. Suppose this is the case. 

Then one gets after some rearrangement of variables 

\e(iS{z) I 

where S{z) C E{u, v) denotes a set of edges e such that e is labeled as G E{v) and = 1. Substituting 
it into the integral Eq. (1381 ). taking into account that 6{e) is a central element and that J D6{e) 0(e) = 1 one 
gets 

Tu.. = i0u,ir' ■ ■ ■ {0u,pr' iOv,b+iy' • • • iOvM.^' 

which coincides with the desired expression Eq. (|37] ). □ 

Corollary 4. Suppose and are matchgates. Then the contracted tensor T^i^v cil^o a matchgate. 

Proof. Since cyclic shifts of indexes map matchgates to matchgates, see Corollary |2]in Section [X2l we can 
assume that the edges of T„ and T„ are already ordered as required in Lemma [5] Represent T^, by their 
canonical generating functions, see Theorem |2l Using Eq. (|38] ) one concludes that Tui,v{0) is a Gaussian 
integral I{A, B) for some matrices A and B, see Eq. ([TTI ). Therefore, T^*^ is a matchgate, see Corollary [3] 
in Section [331 □ 

Remark: Given the canonical generating functions for Tu and T^, the canonical generating function for 
the contracted tensor T^^^ can be obtained straightforwardly using Eq. ([38]) and computing the resulting 
Gaussian integral I{A, B) using Eq. ([141 ). The details can be found in Appendix A. 
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Figure 6: An open tensor network with 7 external edges equipped with a Kasteleyn orientation. 



4.3 Contraction of a planar subgraph in one shot 

Suppose a planar connected graph G = {V, E) is a part of a larger non-planar tensor network such that G is 
connected to the rest of the network by a subset of external edges Eext C E. The remaining internal edges 
Eint = E\Ef.xt are the edges that can be contracted "locally" without touching the rest of the network. By 
abuse of definitions, we shall assume that the external edges have only one endpoint (the other endpoint 
belongs to the rest of the network) which belongs to the outer face of G, see Fig. 6. For convenience let us 
also assume that the graph G is embedded into a disk such that the external edges stick out from the disk as 
shown on Fig. 6. A network that consists of such a graph G = {V, E) and a collection of tensors {Tu}u<^v 
will be referred to as an open tensor network. Throughout this section we shall consider only open tensor 
networks in which every tensor is a matchgate. Contraction of an open tensor network amounts to finding 
a tensor Ty of rank \Ef.xt\ obtained by contracting all internal edges of G. It follows from Corollary |4l 
Section 14.21 that Ty is a matchgate. The goal of the present section is to represent the generating function 
for the contracted tensor Ty as a convolution integral similar to Eq. (1381 ) where the integration is taken over 
all internal edges. 

An alternative strategy for computing Ty is to apply the edge contraction described in the previous section 
sequentially until all internal edges of G are contracted. Although it yields a polynomial-time algorithm this 
strategy is not very robust. An obvious drawback is that every edge contraction involves computing the 
Gaussian integral Eq. (fT4l) which requires a matrix inversion. Contracting sequentially 0{n) edges would 
require 0{n) nested matrix inversions which may be difficult or impossible to do if the matrix elements 
are specified with a finite precision. In order to reduce the number of nested matrix inversions one could 
organize the edge contractions into a sequence of rounds such that each round involves contractions of 
pairwise disjoint edges. The contractions involved in every round can be performed in parallel. The number 
of the rounds can be made O(logn) using the techniques developed by Fiirer and Raghavachari 1221 . We 
shall not pursue this strategy though because the approach described below allows one to compute Ty using 
only one matrix inversion. 

The main result of this section is the following theorem. 

Theorem 4. Consider an open matchgate tensor network on a planar graph G = (y, E) with n vertices and 
m external edges. Assume that the tensors Ti, . . . , T„ are specified by their canonical generating function, 



Then the tensor Ty obtained by contracting all internal edges ofG can be represented as a Gaussian integral 





(39) 
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Here A, B are matrices of size kxk and k x mfor some k = 0{{n + m)^). Matrix elements of A and B are 
linear functionals of Ai, . . . , An and Bi, . . . , Bn- One can compute A and B in time 0{k). Furthermore, 
ifG has bounded vertex degree then the same statement holds for k = 0{m + n). 

Before going into technical details let us explain what is the main difficulty in representing the contracted 
tensor Ty by a single Gaussian integral. The point is that the convolution formula Eq. (l38l) holds only if 
the edges incident to the vertices n, v are ordered in a consistent way as shown on Fig. 4. If the orderings 
are not consistent, an extra sign may appear while commuting the variables living on the contracted edges 
towards each other. Assume one wants to contract the combined vertex u-kv with some third vertex w. If the 
ordering of edges at the combined vertex n ★ v is not consistent with the ordering at w, one has to perform 
a cyclic shift of indexes in the tensor Tui,v and/or before one can directly apply the formula Eq. (l38l) to 
Tui,v and T^. Therefore, in general one can not represent the tensor Tui^v-kw obtained by contracting u, v, w 
as a single Gaussian integral. 

In order to avoid the problem with inconsistent edge orderings we shall contract an open matchgate tensor 
network in two stages. At the first stage one simulates each tensor by a matching sum of some planar 
graph as explained in Section 13.31 It yields an open tensor network in which every tensor has a linear 
generating function (since every vertex must have exactly one incident edge). At the second stage one 
represents the contraction of such a network by a single convolution integral analogous to Eq. (l38l) . The 
problem with inconsistent edge ordering will be addressed by choosing a proper orientation on every edge 
(which affects the definition of monomials 9{e) in Eq. (l38l)). One can regard this approach as a generalization 
of the original Kasteleyn's method [8] to the case of a matching sum with "boundary conditions". 

Definition 3. A tensor T is called linear if it has a linear generating function, T = ^^1^=1 "^a Qa- 

Clearly, any linear tensor T can be mapped to T{9) = 0i by a linear change of variables. Lemma [l] 
implies that T{6) = 6*1 is a matchgate. Therefore any linear tensor is a matchgate, see Lemma[2l 



Definition 4, Orientation of a graph G = iy, E) is an antisymmetric matrix A of size \V\ x \V\ such that 

A^i 



'■u,v 



±1 if iu,v)(£E, 
otherwise. 



An edge (u, v) E is oriented from u to v iff A^^v = 1- 

Recall that we represent each tensor T„ by a generating function Tu{9) that depends on Grassmann 
variables {Ou^i, • • • , Ou,d{u)) associated with the edges incident to u, see Eq. (l36l ). Given an orientation A of 
the graph G and an edge e = {u,v) e E with the labels e" G E{u) and G E{v), define 



9{e) 



/ D9{e) = Au,v / d9,,k / d9u,j, and / De{e) = \{ De{e 



(40) 

Note that 9{e) and / D9{e) are symmetric under the transposition of u and v. 

Lemma 6. Let Ty be a tensor obtained by contraction of an open tensor network on a graph G = (y, E). 
Assume that all tensors in the network are linear Then there exists an orientation A and an ordering of the 
vertices V = {fi, f2, . . . ,Vn} such that 



I D9{e)n,n,---n„ exp [ ^ ^(e) I . 



Tv = / D9{e) • • • exp I > ^ ^(e) | . (41) 

The orientation and the ordering can be found in time 0{n). 
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Remark 1: The generating function of Ty is defined for the ordering of the external edges in which they 
appear as one circumnavigates the boundary of the disk anticlockwise. The order of variables in Ty corre- 
sponds to the counterclockwise order of the external edges. 

Proof. Without loss of generality G is a 2-connected graplH. Then the boundary of the outer face of G is a 
closed loop without self-intersections. Let us denote it Tout- Mark some vertex in Tout that has at least one 
incident external edge (if there are no external edges, mark an arbitrary vertex). Let Tout = {1> 2, . . . , m} be 
an ordered list of all vertices on the outer face of G corresponding to circumnavigating Tout anticlockwise 
starting from the marked vertex. Extend the ordering of vertices to the rest of V in an arbitrary way, so that 
V = {1,2, . . . ,n} and the first m vertices belong to Tout- 

Definition 5. Let G be a planar graph with the vertices ordered as described above. A Kasteleyn orientation 
(KO) ofG is an orientation A such that 

(1) The number ofc.c.w. oriented edges in the boundary of any face of G is odd (except for the outer face). 

(2) ^1,2 = ^2,3 = • • • = Am-l,m = 1- 

Remark: The standard definition of a KO requires that (1) holds for all faces of G including the outer face 
and does not require (2), see for example llT3l . By abuse of definitions we shall apply the term KO to 
orientations satisfying (1),(2). The standard definition is not suitable for our purposes because G may have 
odd number of vertices while the standard KO exists only on graphs with even number of vertices. The 
condition (2) is needed to ensure consistency between different "boundary conditions". Example of a KO is 
shown on Fig. 6. 

Proposition 3. Any planar graph has a KO. It can be found in a linear time. 

We postpone the proof of the proposition until the end of the section. Let us choose the orientation A 
in Eq. (l40l ) as a KO of the graph obtained from G by removing all external edges. Let us verify that the 
contracted tensor Ty can be written as in Eq. (|4T]) . 

Indeed, let S C E^xt be a subset of external edges such that any vertex in {1, ... , m} has at most one 
incident edge from S. (Below we shall consider only such sets S without explicitly mentioning it.) Let dS 
be a set of vertices that have an incident edge from S (clearly all such vertices belong to the outer face). For 
any S as above and any 95-imperfect matching M € M-{G, dS) define a subset of Grassmann variables 

Q{S, M) = {(n, j) : ueV, and e| G 5 U M}. 

In other words, (u, j) € Q{S, M) iff Ouj is a Grassmann variable that live on some edge of S" U M. Note 
that there are two Grassmann variables living on any internal edge and one variable living on any external 
edge. Thus for any S and M the set Q.{S, M) contains n variables. Define a normally ordered monomial 

n Gu,3 (42) 

{u,j)&n{s,M) 

as a product of all variables in Q{S, M) ordered according to 

(01^1, . . . ,6*1,^(1), 612,1, • • • , ^2,d(2); • • • ) ^n,l, • • • ,On,d(n))- (43) 

Define also M-ordered monomial 

n n ^(^)' (44) 

(«,i):e^eS' eeM 

^If G has a cut- vertex it one can always add an extra edge to some pair of nearest neighbors of u in order to make G 2-connected. 
The new edge must be assigned a zero weight in the two tensors it belongs to. Since the new edge does not contribute to Tv it can 
be safely removed at the end of the analysis. 
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where the order in the first product must agree with the chosen ordering of edges in Egxt, see Fig. 6. Clearly 
the two products Eqs. (|42I44I) coincide up to a sign that we shall denote sgn(M). In order to prove Lemma[6] 
it suffices to show that 

sgn(M) = 1 for all 95-imperfect matchings M, for all S C E'g^.^. (45) 
Indeed, denoting T„ = Yl'j=i ^u,j one can rewrite Eq. (|4TI ) as 

Tv = Yi E / ^^(^) n <^-jn^(^) 

SCEe:ct MeMiG,dS) '^'"^^'^"^ {u,j)en(S,M) e^M 

= E n E 'Sn{M) H wl (46) 

Assuming sgn(M) = 1 one can identify the sum over M € M.{G, dS) with the component of the contracted 
tensor Ty in which the subset S of external edges carries index 1. 

Note that for any S C E^xt and any 95-imperfect matching M each vertex u £ V contributes exactly 
one variable to Q,{S, M). Indeed, at every vertex u £ V there is either one incident edge from M or one 
incident external edge. All other edges incident to u and the variables living on these edges can be ignored 
as far as computation of sgn(M) is concerned. Therefore one can compute the sign sgn(M) by introducing 
auxiliary Grassmann variables rj = (ryi, . . . , r?„) associated with vertices of G and comparing the normal 
ordering of ( the one in which the indexes increase from the left to the right) with the M -ordering of r], 
namely 

II II "^^^^ ^ sgn(M) r/i??2 • • • ??n, where 77(e) = A^^v VuVv if e = {u, v). 

u&dS eeM 

Here the ordering in the first product is normal while the ordering in the second product may be arbitrary 
since 7]{e) is a central element. Consider any subsets S, S' C E^xt- Given any (9S'-imperfect matching M 
and c?5"-imperfect matching M' define a relative sign 

sgn(M, M') = sgn(M) sgn(M'), (47) 

such that 

n n r]{e) = sgn{M,M') Vu v{e). (48) 

u€dS e&M u£dS' eGM' 

In order to compute sgn(M, M') consider the symmetric difference M M'. It consists of a disjoint union 
of even-length cycles Ci, . . . , Cp and open paths Fi, . . . , Fg such that every path Tj has both its endpoints 
in the symmetric difference dS © dS'. Given a path Tj with endpoints s, t G dS dS', s < i let us orient 
Tj from s to t. Now one can compute the relative sign as follows. 

Proposition 4. Consider any subsets S, S' C E^xt- Let Ci, . . . ,Cp and Fi, . . . , F^ be the cycles and the 
paths formed by M © M' for some dS -imperfect matching M and some dS' -imperfect matching M'. For a 
path Tj connecting vertices s,t € dS dS' on the outer face such that s < t let uj{Tj) = 1 if the interval 
(s, t) contains odd number of vertices from dS and io{Tj) = if this number is even. Then 

sgn(M,M') = (-l)P n$(C,) n(-ir(r'=)$(Ffc), (49) 
i=i k=i 

where 

HCj)= n H^k)= n 
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Remark 1: The definition of u}{Tj) is symmetric under exchange of S and S' . Indeed, the overall number 
of vertices from dS © dS' contained in the interval (s, t) is even since these vertices are pairwise connected 
by r's. The remaining vertices of (s, t) either belong to both sets S, S' or belong to neither of them. 
Remark 2: The product t;)Grfe -^u.v gives the parity of the number of edges in whose orientation 
determined by A disagrees with the chosen orientation of T^. The product ^{Cj) does not depend on how 
one chooses orientation of Cj since every cycle Cj has even length. 

Proof. Indeed, one can easily check that for every cycle Cj one has 

n r?(e) = -$(C,) n ^(^)- (50) 

Therefore changing the M-ordering to the M'-ordering in a cycle Cj contributes a factor — <I>(Cj) to the 
relative sign sgn(M, M'). Consider now a path Vj connecting vertices s, t € dS © dS' where s <t. Let us 
argue that changing the M-ordering to the M'-ordering on the path Tj contributes a factor (— l)^^"^) (^{Tj) 
to the relative sign sgn(M, M'). Indeed, one can easily check the following identities: 

: iism n ^(e) = ^(ri) n ^(^)' 

s,t £ S' : the same as above up to M ^ M', 
seS,teS' : 7?, II 7?(e) = cD(r,)r/t J] r?(e), 

s £ S' ,t G S : the same as above up to M <-> M'. 

Consider as example the case s,t G S. Bringing the variables rjs and rjt together in the monomial Hweas Vu 
introduces an extra sign (— 1)'^'^^^^ Taking into account that r]{e) are central elements and using the first 
identity above one concludes that 

n vie) = {-ir^^^^ HT,) n ^« n ^(^)- 

u&dS eeVjCiM u(^dS\{s,t} eeV^nM' 

Other three cases can be considered analogously using Remark 1 above. Combing it with Eq. (l50l) one 
arrives to Eq. ([491). □ 



Let us proceed with the proof of Lemma [6] The first condition in the definition of KO implied that 
$(Cj) = —1 for all cycles Cj. Indeed, consider any particular cycle Cj and let Nq, Ni,N2 be the number of 
vertices, edges, and faces in the subgraph bounded by Cj. The Euler formula implies that NQ-'rN2 — Ni = 1. 
Denote also A^|"* the number of internal edges, i.e., edges having at least one endpoint in the interior of Cj. 
Since Cj has even length, iV|"* has the same parity as A'^i. Furthermore, since all vertices of the subgraph 
bounded by Cj are paired by M (and by M'), Nq is even. Since $(Cj) can be regarded as a parity of c.c.w. 
oriented edges in Cj and each internal edge is c.c.w. oriented with respect to one of the adjacent faces the 
property (1) of KO yields 

Therefore Proposition |4] implies 

1 

sgn(M,M') = Wi-ir^^^^^Tk). (52) 

k=l 



^This is the well-known property of a Kasteleyn orientation which we prove below for the sake of completeness. 



24 



Let us now show that 

(_l)-(r.)$(r^) = 1 (53) 

for all paths Ffc. Indeed, let s, t € 5 © S' be the starting and the ending vertices of Fk- Consider a path 
obtained by passing from t to s along the boundary of the outer face Tout in the clockwise direction. Let 
A'o, Ni,N2 be the number of vertices, edges, and faces in the subgraph bounded by a cycle U F^. Denote 
also Nf^*' the number of edges that have at least one endpoint in the interior of Tjt U F^. The Euler formula 
implies that A'^o + — -^i = 1- Note that ^(F^) can be regarded as the parity of the number of edges in 
Fjt whose orientation determined by A corresponds to c.c.w. orientation of the cycle F^ U F^,. Repeating the 
arguments leading to Eq. (|5TI) and noting that all edges of the cycle F^ U F^ belonging to F^ are oriented 
c.c.w. one gets 

$(Ffe) = (-i)ini+^2+7vr* ^ ^_^^ir,\+N2+Ni ^ ^_-^^|rfe|+iVo+i_ (54) 

Here \Tk \ and |F^| are the numbers of edges in the two paths. Consider three possibility: 

Case 1: s,t G dS. Then |Ffc| is odd and thus $(Ffc) = All Nq vertices of the graph bounded by 

Ffc U F^ are paired by the matching M except for s,t and those belonging to dS and lying on the interval 

(s, t). Therefore the parity of A'^o coincides with w(Ffc) and we arrive to Eq. (1531 ). 

Case 2: s, t G dS'. The same as Case 1 (see Remark 1 after Proposition HJl. 

Case 3: s £ dS,t e dS' (or vice verse).Then [F^l is even and thus ^{Tk) = (-1)^°+^ All Nq vertices 
of the graph bounded by F^ U F^ are paired by the matching M except for s (or except for t) and those 
belonging to dS and lying on the interval (s, t). Therefore the parity of Nq coincides with uj{Tk) + 1 and 
we arrive to Eq. (1531 ). 

Combining Eqs. (15 11531 ) and Proposition |4] we conclude that sgn(M, M') = 1 for all M and M'. Thus 
either sgn(M) = 1 for all M or sgn(M) = —1 for all M. One can always exclude the latter possibility by 
applying a gauge transformation to the orientation A. A gauge transformation at a vertex u £ V reverses 
orientation of all edges incident to u. Let us say that a vertex n G 1/ is internal if does not belong to the 
outer face of G. Clearly a gauge transformation at any internal vertex u maps a KO to a KO and flips the 
sign sgn(Af ) for all M. Thus it suffices to consider the case when G does not have internal vertices (i.e. 
G is an outerplanar graph). If m = n is even, a matching M = {(1, 2), (3, 4), . . . , (m — 1, m)} has sign 
sgn(M) = 1 due to property (1) of a KO and thus all matchings have sign +1. If m = n is odd one can 
apply the same argument using a matching M = {(2, 3), (4, 5), . . . , (m — 1, m)} (recall that the vertex 1 
has at least one external edge and thus it can be omitted in M). □ 

Proof of Theorem^ Let rig be the number of internal edges in the graph G, so that \E\ = m + Ue- Since 
G is a planar graph, = 0{n), see for example [24J, and thus \E\ = 0{n + m). Denote degree of a 
vertex u G F by d{u) (it includes both internal and external edges). Applying Theorem[3]one can simulate 
the tensor at any vertex u G F by a matching sum of some planar graph Gu with 0{d{uf') vertices. 
Combining the graphs Gu together one gets an open tensor network G' = {V ,E') in which all tensors are 
linear. The network G' has m external edges. The number of vertices n' in the network G' can be bounded 
as n' = Y^^^y 0{d{uf) = 0{{J2uev diu)^) = 0{\E\'^) = 0{{m + nf). If G has bounded degree one 
gets n' = X^^gy 0{d{uf') = 0{n). Thus in both cases n' = 0{k), where k is defined in the statement of 
the theorem. It follows from Theorem [3] that the edge weights in the matching sums are linear functions of 
the matrix elements of yli, . . . , A„ and . . . , Let n'^ be the number of internal edges in G' . Since G' 
is a planar graph, n'^ = O(n') = 0{k). Thus the total number of edges in G' is |E"| = n'^ + m = 0{k). 
Invoking Lemma[6]we need to introduce a pair of Grassmann variables for every internal edge of G' and one 
variable for every external edge. Thus the total number of Grassmann variables is 0{k). It determines the 
number of variables in the vector 9 in Eq. ( [39l ). Representing linear tensors Tj as Gaussian integrals, namely 



Tj = J dfiexp {fJ.Tj), 
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Figure 7: Left: a tensor network with a single vertex embedded into a torus. Right: the pairing graph P. 



one can combine the multiple integrals in Eq. (|4TI ) into a single Gaussian integral Eq. ( [39l ) with the matrix 
A having a dimension 0{k) x 0{k) and B having a dimension 0{k) x m. Thus vl and B have the desired 



Proof of Proposition\3\ Let G = (F, E) be a planar graph with n vertices such that the outer face of G is a 
simple loop. An orientation A satisfying (1) can be constructed using the algorithm of [TT]. For the sake of 
completeness we outline it below. Let G* = {V* ,E) be the dual graph such that each face of G contributes 
one vertex to G* (including the outer face). Let T be a spanning tree of G* such that the root of T is the 
outer face of G. One can find T in time 0(|y| + \E\) = 0{n) since for planar graphs \E\ = 0(|y|). 
Assign an arbitrary orientation to those edges of G that do not belong to T. By moving from the leaves of 
T to the root assign the orientation to all edges of T. Note that for every vertex u of T which is not the 
root the orientation of an edge e connecting u to its ancestor is uniquely determined by (1). We obtained an 
orientation of all edges of G satisfying (1). 

In order to satisfy (2) one can apply a series of gauge transformations. A gauge transformation at a 
vertex u ^ V reverses orientation of all edges incident to u. Clearly it preserves the property (1). Applying 
if necessary a gauge transformation at the vertices {1, 2, . . . , m — 1} one can satisfy (2). □ 

4.4 Contraction of matchgate networks with a single vertex 

In this section we explain how to contract a matchgate tensor network T that consists of a single vertex 
u with m self-loops embedded into a surface S of genus g without self-intersections. Example of such a 
network with m = 3 and g = lis shown on Fig. 7. Let T be a tensor of rank 2m associated with u. Clearly 
the contraction value c(T) depends only on the pairing pattern indicating what indexes of T are contracted 
with each other. It will be convenient to represent the pairing pattern by a pairing graph P = {V, E) with a 
set of vertices V = {1,2, . . . , 2m} such that a pair of vertices (a, b) is connected by an edge iff the indexes 
o, 6 of the tensor T are contracted with each other (connected by a self-loop). By definition P consists of m 
disjoint edges. Let us embed P into a disk such that all the vertices of P lie on the boundary of the disk and 
their order corresponds to circumnavigating the boundary anticlockwise. The edges of P are represented 
by arcs lying inside the disk, see Fig. 7. One can always draw the arcs such that there are only pairwise 
intersection points. 

Introduce an auxiliary tensor R of rank 2m such that 



properties. 



□ 




1 if Xa = Xb for all (a, b) G E, 
if Xa Xf) for some (a, b) G E. 



The contraction value of T can be represented as 





(55) 
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Let 9 = (6*1, ... , 92m) and -q = {rji, . . . , ri2m) be Grassmann variables and T{9), R{i]) be the generating 
functions of T and R. 

Proposition 5. Let e(T) = 0, 1/or even and odd tensors T respectively . Then 

c{r) = i<^^ j D9 j DriT{9)R{i]) exp(i0^7/). (56) 

Proof. A non-zero contribution to the integral comes from the terms in which T{9) contributes monomial 
T{x) 9{x) and R{'q) contributes monomial R{x) r]{x) for some x G {0, 1}^™. A simple algebra shows that 
for any x G {0, l}^*" one has the following identity 

9{x)v{x) J] iMa = rN(-l)N(l^l-i)/2e(l2-)r/(l2-), 

a : Xa =0 

where \x\ is the Hamming weight of x. Taking into account that T(x) = unless |x| has parity e(T) one 
gets 

(_l)kl{l^-|-l)/2 = ^-€(T)_ 

Since J D9 J D^9{1^'^) r/(l2™) = 1, one gets Eq. □ 

In general R is not a matchgate tensor because the chosen planar embedding of the pairing graph may 
have edge crossing points. For example, assume that P has 4 vertices {1,2,3,4} and two edges (1,3), 
(2, 4) (which can be realized on a torus). Then the non-zero components of R are i?(0000) = i?(1010) = 
i?(0101) = i2(llll) = 1. Substituting them into the matchgate identities Eq. (|2]) for even rank;-4 tensors 
one concludes that R is not a matchgate. 

Let us order the edges of P in an arbitrary way, say, E = {ei, 62, ... , Cm}- For any edges Cp, Cq & E 
let Np^q be the the number of self-intersections of ep, tq in the planar embedding shown on Fig. 7. Since 
we assumed that all intersections are pairwise, Np^q takes only values 0, 1, i.e., is a symmetric binary 
matrix. Let us also agree that Np.p = 0. We shall see later that the tensor R can be represented as a linear 
combination of 2*" matchgate tensors, where r is a binary rank of the matrix A^. It is crucial that the rank of 
A^ can be bounded by the genus g of the surface S. 

Lemma 7. The matrix N has binary rank at most 2g. 

Proof. Let us cut a small disk D centered at the vertex u from the surface S, embed the pairing graph P 
into the disk D as shown on Fig. 7 and glue the disk back to the surface S. Thus given any self-loop a 
connecting indexes a and b of the tensor T, a small section of a lying inside D is replaced by an edge 
e = {a,b) £ E of the pairing graph. We get a family of m closed loops embedded into S. The loops may 
have pairwise intersection points inside the disk D. Let Op be a loop that contains an edge Cp G E. To every 
loop Up one can assign its homological class [op] G //i(S,Z2) in the first homological group of S with 
binary coefficients. Since all intersection points between the loops are contained in the disk D, we get 

Np^q = uj{[ap], [aq]), 

where uj : Z2) x Z2) ^ {0, 1} is the intersection form. It is well known that the intersection 

form defined on a surface S of genus g has rank 2g. Therefore, A^ has rank at most 2g. □ 

Given any edge e £ E, let Z(e),r(e) G F be the two endpoints of e such that /(e) < r(e). Denote 
'q{e) = rji(^e) Vr{e)- The generating function for the tensor R can be written as 

Riri)= Yl (-l)^'^'^''n'?('^)' ^here vie) = m(e)Vriey (57) 
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Here we identified a binary string y G {0, 1}'" with the subset of edges ea £ E such that Ua = 1- Indeed, 
for any x G {0, l}^™ such that R{x) = 1 one has to regroup the factors in r]{x) to bring together variables 
corresponding to the same edge. It yields an extra minus sign for every pair of intersecting edges in y. Since 
every pair of edges Cq, Cb contributes a sign (^—I'j^a.byaVb^ arrive to Eq. (|57]) . 

Consider binary Fourier transform of the function (—1)2^ ^ ^, 

/W = ^ E (-1)^^"^^+^-^ z€{0,ir. (58) 

S/G{0,1}™ 

Clearly f{z) = unless z G Kev{N)-^, where Ker(iV) = {y G {0, 1}™ : iVy = 0} is the zero subspace 
of N. If N has rank r, the zero subspace of N has dimension m — r and thus Ker(A^)-'- has dimension r. 
Let us order all the vectors of Ker(A^)^ in an arbitrary way 

Ker(iV)^ = 

Applying the reverse Fourier transform one gets 

{-i)W^y = Y^f{z-){-i)y^\ (59) 

a=X 



By Lemma|7]the number of terms in the sum above is bounded by 2^^. Substituting Eq. (|59l ) into Eq. ( [57] ) 
we arrive to 

R{ri) = E /(^'^) ( E(-l)^'"^^ ^(^) ) ' (60) 

a=l VeG-E / 

where {z"')e is the component of the vector corresponding to an edge e. It shows that R is indeed a linear 
combination of 2*" matchgate tensors with r < 2g. 

In order to get an explicit formula for the contraction value Eq. ( [55] ) let us introduce an auxihary 2m x 2m 
matrix 

+1 if j = l{e), k = r(e) for some e £ E, 
Aj^k = \ —1 if J = k = /(e) for some e G -E 
otherwise 



(_l){2")e if j = foj- some e e E, 
1 otherwise. 



Introduce also auxiliary diagonal 2m x 2m matrices D"^, a = I, ... ,2^ such that 

Then Eq. ([60] ) can be rewritten as 

R{^) = E •^(^") A I)'^ r/) . (61) 

Theorem [2] implies that T can be described by a generating function 

r(e) = C exp Q 0^ F 0^ Dfi exp (/x^ (7 0) , 

where F and G have size 2m x 2m and k x 2m, for some even integer < < 2m. Using Eq. ([56] ) one 
can express the contraction value c(T) as a linear combination of 2^' Gaussian integrals 

c(T) = C ^/(z'^) J DeDriDfiexp(^^e'^F9+^T]^D''AD''ri + fi^G9 + ie^r]j. (62) 
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Introducing a (4m + k) x (4m + k) matrix 



one finally gets 



F il -G^ 

M" = -il -D^AD" 
GOO 



2'' 

c(T) = C J]/(z") Pf(M"). 

a=l 



(63) 



Computing Pf {M°-) requires time O(m^). Lemma|7]implies that the number of terms in the sum is at most 
2^3 _ Finally, as we show below one can compute f{z"') in time O(m^). Thus c(T) can be computed in time 
0(m3)229. 



Proposition 6. The function f(z) in Eq. d58D can be represented as 

1 



for some matrix M computable in time 0{m^). 



2r/2 



(-1) 



hz^Mz 



(64) 



Proof. Using Gaussian elimination any symmetric binary matrix N with zero diagonal can be represented 
as = N U, where f7 is a binary invertible matrix and iV is a block diagonal matrix with 2x2 blocks, 

r/2 



"-e(?i). 

i=i ^ ^ 



In particular, the rank of is always even. The matrix U can be found in time O(m^). Performing a change 
of variable y Uy in Eq. (1581 ) one gets 



(65) 



ye{0,l}'" 

It follows that f{z) = unless Zr+i = . . . = Zm = 0. Using an identity 

(_l)^i'^2 _ 1 ^_iYi-y2+yi-xi+y2-x2 

^ S/l,?/2=0,l 



one can rewrite Eq. (1651) as 



2r/2 ^ ^ 2''/2 

We get the desired expression Eq. ([641) with M = U^^ N (U^^)'^ 



□ 



4.5 The main theorem 

Theorem [T] can be obtained straightforwardly from Theorem |4] and the contraction algorithm for a network 
with a single vertex, see Section 14.41 Indeed, let M be a planar cut of G with m edges and Gm be a 
subgraph obtained from G by removing all edges of M. By definition Gm is contained in some region D 
with topology of a disk. Without loss of generahty D contains no edges from M (otherwise one can remove 
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these edges from M getting a planar cut with a smaller number of edges). Thus one can regard Gm as an 
open tensor network with 2m external edges. Since Gm contains all vertices of G, the network obtained by 
contraction of Gm consists of a single vertex and m self-loops. As explained in the previous section, one 
can compute the contraction value of such a network in time O(m^) 2^^. 

In order to contract Gm one has to compute the Gaussian integral Eq. (|39l ). Theorem H] guarantees that 
this integral involves matrices of size k, where k = 0{{n + mY) or k = 0{n + m) depending on whether 
the graph G has bounded vertex degree. As explained in Section 1231 the Gaussian integral with matrices of 
size k can be computed in time 0{k^). Combining the two parts together one gets Theorem[T] 
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Appendix A 

Suppose Tu and T„ are matchgate tensors specified by their canonical generating functions as in Eq. ([TT] ). 
that is 



Toe = Coc exp Q e'l Aa 6a^ j exp {nl Ba 9a) , where a 



U, V. 



Here 6u = • • • , and 6^ = {O^^i, . . . , are the two sets of Grassmann variables associated 

with the vertices u and v. Denote also e{T) the parity of a matchgate tensor T, that is, e{T) = (e(T) = 1) 
for even (odd) tensor T. In the remainder of this section we explain how to express the canonical generating 
function for the contracted tensor Tui,v, see Eqs. (I37I38I) . in terms of the matrices Aa, Ba- 
Applying Eq. (1381 ) one gets 

(66) 



GuG, / De{e) / Dfiu / Dfi 

V 6Xp [f{9ui Svi 

Je£E{u,v) J J 

f{0u,O^,^iu,^^v) = Yl l€AaOa + fllBaea+ Yl 



where 

f(f)...f),..ii„..ii„.) = V 

a=u,v e£E{u,v) 

Let us split the vectors of Grassmann variables 6'„, 9y into external and internal parts. 



0u — {Ou,O\^) and 9^ — {91,91), 

such that all internal variables are integrated out in T^^^,. Then one can rewrite the expression Eq. ( [66l ) as 
a product of a Gaussian exponent and the standard Gaussian integral I{K, L), see Eqs. (I13I14I ). for some 
matrices K, L defined below, 

Tu.v{r) = C„C,(-l)^+<^")^(^") exp Q^^^^) / ^^exp Q r/^ r? + r?^ L . (67) 

Here we introduced auxiliary vectors of Grassmann variables r = (0^,0^), r/ = (0^^, 0* , /i^, /i^,). The 
matrices H, K, L above will be defined using a partition of matrices Aa, Ba into "internal" and "external" 
blocks as follows: 



Aee Aei 
Aie Ail 



Ay 



All Aie 
Aei Aee 



Bu=[ Bl -Bu ] , By=[Bl Bl] . 
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Introduce also a square matrix I that has ones on the diagonal perpendicular to the main diagonal and zeroes 
everywhere else. Then the matrices H, K, L in Eq. (l67l) are defined as 

■ 

Ai^ 
Bl • 
Bl 

Finally, the extra sign in Eq. (l67l) takes into account the difference between the order of integrations in 
Eqs. (I66I67I) . Summarizing, Eq. (|67] ) together with the Gaussian integration formulas Eqs. (I13I14I) allow one 
to write down the canonical generating function for the contracted tensor T^i^v 
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